------------------------------------------------------------------------------------------------------------------------------------------
      name:  <unnamed>
       log:  C:\Users\mdnieman\Desktop\SS Rep\SS_SMPO_epa.log
  log type:  text
 opened on:  29 May 2018, 16:07:52

. /*******************************/
> /*  Mark David Nieman  */
> /*  1/21/2016          */
> /*  Application: EPA data */
> /*  Stata 12.1        */
> /*******************************/
> 
> #delimit;
delimiter now ;
. clear;

. clear matrix;

. clear mata;

. program drop _all;

. set seed 1999;

. /* strategic logit with partial observability */
>         #delimit;
delimiter now ;
.         program define slpo_lf, rclass;
  1.                 args lnf A_L A_Rr B_Rr;
  2.                 tempvar pB_r pA_R ;
  3.                         quietly gen double `pB_r' = (exp(`B_Rr'))/(1+exp(`B_Rr'));
  4.                         quietly gen double `pA_R' = (exp(`pB_r'*`A_Rr'))/(exp(`A_L')+exp(`pB_r'*`A_Rr'));
  5.                                         quietly replace `lnf' = ln((1-`pA_R') + (`pA_R')*(1-`pB_r')) if $ML_y1==0;
  6.                 quietly replace `lnf' = ln((`pA_R')*(`pB_r')) if $ML_y1==1;
  7.         end;

.         /*load data*/
> use CAA_Dataset1.dta, clear;

. /* set=1 if pp or hospital*/
> 
> sort afsid year;

. gen id=_n if afsid!=afsid[_n-1];
(650,472 missing values generated)

. replace id=id[_n-1] if afsid==afsid[_n-1];
(650,472 real changes made)

. drop if id==.;
(0 observations deleted)

. drop if afsid==afsid[_n-1] & year==year[_n-1];
(0 observations deleted)

. isid id year;

. tsset id year;
       panel variable:  id (strongly balanced)
        time variable:  year, 1999 to 2011
                delta:  1 unit

. /* logit analysis */
> logit nc_year major  public PercentBlack PercentHisp PercentPoverty MedHsInc PercentUnemploy L.nc_year;

Iteration 0:   log likelihood = -196872.24  
Iteration 1:   log likelihood =  -160340.6  
Iteration 2:   log likelihood = -113820.16  
Iteration 3:   log likelihood = -111873.42  
Iteration 4:   log likelihood = -111859.61  
Iteration 5:   log likelihood =  -111859.6  

Logistic regression                             Number of obs     =    650,472
                                                LR chi2(8)        =  170025.27
                                                Prob > chi2       =     0.0000
Log likelihood =  -111859.6                     Pseudo R2         =     0.4318

---------------------------------------------------------------------------------
        nc_year |      Coef.   Std. Err.      z    P>|z|     [95% Conf. Interval]
----------------+----------------------------------------------------------------
          major |   .8070315   .0120971    66.71   0.000     .7833217    .8307413
         public |   .0114424   .0265084     0.43   0.666     -.040513    .0633979
   PercentBlack |    .000432   .0005458     0.79   0.429    -.0006378    .0015018
    PercentHisp |   .0043885   .0005294     8.29   0.000     .0033508    .0054262
 PercentPoverty |  -.0195227   .0022273    -8.77   0.000    -.0238882   -.0151572
       MedHsInc |  -6.67e-07   8.33e-07    -0.80   0.423    -2.30e-06    9.65e-07
PercentUnemploy |  -.0027828   .0025835    -1.08   0.281    -.0078463    .0022806
                |
        nc_year |
            L1. |   4.371479   .0126605   345.28   0.000     4.346664    4.396293
                |
          _cons |  -3.423012   .0592595   -57.76   0.000    -3.539159   -3.306866
---------------------------------------------------------------------------------

. matrix ll_l = e(ll);

.         gen sample = 1 if e(sample)==1;
(54,206 missing values generated)

. sum nc_year if sample==1;

    Variable |        Obs        Mean    Std. Dev.       Min        Max
-------------+---------------------------------------------------------
     nc_year |    650,472    .0900531     .286258          0          1

. tab nc_year if sample==1;

      (max) |
    nc_year |      Freq.     Percent        Cum.
------------+-----------------------------------
          0 |    591,895       90.99       90.99
          1 |     58,577        9.01      100.00
------------+-----------------------------------
      Total |    650,472      100.00

. sutex nc_year major  public PercentBlack PercentHisp PercentPoverty MedHsInc PercentUnemploy if sample==1, minmax nobs;
%------- Begin LaTeX code -------%

\begin{table}[htbp]\centering \caption{Summary statistics \label{sumstat}}
\begin{tabular}{l c c c c c}\hline\hline
\multicolumn{1}{c}{\textbf{Variable}} & \textbf{Mean}
 & \textbf{Std. Dev.}& \textbf{Min.} &  \textbf{Max.} & \textbf{N}\\ \hline
nc\_year & 0.09 & 0.286 & 0 & 1 & 650472\\
major & 0.305 & 0.46 & 0 & 1 & 650472\\
public & 0.05 & 0.217 & 0 & 1 & 650472\\
PercentBlack & 11.905 & 14.062 & 0 & 86.489 & 650472\\
PercentHisp & 9.601 & 12.136 & 0.094 & 97.539 & 650472\\
PercentPoverty & 14.85 & 5.509 & 2.117 & 50.888 & 650472\\
MedHsInc & 45055.472 & 11782.858 & 16271 & 119075 & 650472\\
PercentUnemploy & 6.117 & 2.612 & 1.1 & 29.7 & 650472\\
\hline
\end{tabular}
\end{table}
%------- End LaTeX code -------%

. egen N = count(sample) if sample==1;
(54206 missing values generated)

. egen Obs = max(N);

. drop N;

.                         /* Model fit */
>                         /* i's likelihood */
>                         gen double lnfj_l = . ;
(704,678 missing values generated)

.                         replace lnfj_l = ln(invlogit(-(_b[major]*major +  _b[public]*public + _b[PercentBlack]*PercentBlack +
>                                                                                    _b[PercentHisp]*PercentHisp + _b[PercentPoverty]*Perc
> entPoverty + _b[MedHsInc]*MedHsInc + 
>                                                                                    _b[PercentUnemploy]*PercentUnemploy + _b[L.nc_year]*L
> .nc_year + _b[_cons]))) if nc_year==0 & sample==1;
(591,895 real changes made)

.                                                 replace lnfj_l = ln(invlogit( (_b[major]*major +  _b[public]*public + _b[PercentBlack]*P
> ercentBlack +
>                                                                                    _b[PercentHisp]*PercentHisp + _b[PercentPoverty]*Perc
> entPoverty + _b[MedHsInc]*MedHsInc + 
>                                                                                    _b[PercentUnemploy]*PercentUnemploy + _b[L.nc_year]*L
> .nc_year + _b[_cons]))) if nc_year==1 & sample==1;
(58,577 real changes made)

.                         gen double pr_vio_l = invlogit(_b[major]*major  + _b[public]*public + _b[PercentBlack]*PercentBlack +
>                                                                                    _b[PercentHisp]*PercentHisp + _b[PercentPoverty]*Perc
> entPoverty + _b[MedHsInc]*MedHsInc + 
>                                                                                    _b[PercentUnemploy]*PercentUnemploy + _b[L.nc_year]*L
> .nc_year + _b[_cons]) if sample==1;
(54,206 missing values generated)

.                                                 gen pred_vio_l = 0 if sample==1;
(54,206 missing values generated)

.                                 replace pred_vio_l =1 if pr_vio_l > .5 & sample==1;
(52,932 real changes made)

.                         gen corr_l = 0 if sample==1;
(54,206 missing values generated)

.                                 replace corr_l = 1 if (pred_vio_l==1 & nc_year==1) | (pred_vio_l==0 & nc_year==0);
(618,091 real changes made)

.                         gen corr_vio_l = .;
(704,678 missing values generated)

.                                 replace corr_vio_l = 0 if (pred_vio_l==0 & nc_year==1);
(19,013 real changes made)

.                                 replace corr_vio_l = 1 if (pred_vio_l==1 & nc_year==1);
(39,564 real changes made)

.                         gen incorr_vio_l =.;
(704,678 missing values generated)

.                                 replace incorr_vio_l = 0 if (pred_vio_l==0 & nc_year==0);
(578,527 real changes made)

.                                 replace incorr_vio_l = 1 if (pred_vio_l==1 & nc_year==0);
(13,368 real changes made)

.                                         sum corr_l corr_vio_l incorr_vio_l;

    Variable |        Obs        Mean    Std. Dev.       Min        Max
-------------+---------------------------------------------------------
      corr_l |    650,472    .9502192    .2174919          0          1
  corr_vio_l |     58,577    .6754187    .4682222          0          1
incorr_vio_l |    591,895    .0225851    .1485767          0          1

.                                                                                                                                         
>            #delimit;
delimiter now ;
. /* get initial conditions */
> logit nc_year public PercentBlack PercentHisp PercentPoverty MedHsInc PercentUnemploy ;

Iteration 0:   log likelihood = -196872.24  
Iteration 1:   log likelihood = -195910.56  
Iteration 2:   log likelihood = -195898.58  
Iteration 3:   log likelihood = -195898.58  

Logistic regression                             Number of obs     =    650,472
                                                LR chi2(6)        =    1947.32
                                                Prob > chi2       =     0.0000
Log likelihood = -195898.58                     Pseudo R2         =     0.0049

---------------------------------------------------------------------------------
        nc_year |      Coef.   Std. Err.      z    P>|z|     [95% Conf. Interval]
----------------+----------------------------------------------------------------
         public |   .1466835   .0188524     7.78   0.000     .1097334    .1836336
   PercentBlack |  -.0024078   .0003934    -6.12   0.000     -.003179   -.0016367
    PercentHisp |   .0065205   .0003742    17.43   0.000     .0057871    .0072539
 PercentPoverty |  -.0089109   .0015773    -5.65   0.000    -.0120024   -.0058194
       MedHsInc |   6.44e-06   5.74e-07    11.21   0.000     5.31e-06    7.56e-06
PercentUnemploy |   .0490912   .0017545    27.98   0.000     .0456526    .0525299
          _cons |  -2.826456   .0412773   -68.47   0.000    -2.907358   -2.745554
---------------------------------------------------------------------------------

.         matrix B2 = e(b);

.         predict xb;
(option pr assumed; Pr(nc_year))
(54206 missing values generated)

.         gen L = exp(xb)/(1+exp(xb));
(54,206 missing values generated)

.         gen Lmajor =-major;

.         gen Lpp = -pp;

.                 gen R_pub = L*public;
(54,206 missing values generated)

.         gen R_pun = L*L.pun_t;
(54,206 missing values generated)

.         gen R_lag = L*L.nc_year;
(54,206 missing values generated)

.         gen R_ins = L*L.insp_t;
(54,206 missing values generated)

.         gen R_black = L*PercentBlack;
(54,206 missing values generated)

.         gen R_hisp = L*PercentHisp;
(54,206 missing values generated)

.         gen R_pov = L*PercentPoverty;
(54,206 missing values generated)

.         gen R_inc = L*MedHsInc;
(54,206 missing values generated)

.         gen R_unemp = L*PercentUnemploy;
(54,206 missing values generated)

.         gen R_c = L;
(54,206 missing values generated)

. #delimit;
delimiter now ;
. logit nc_year Lmajor R_black R_hisp R_pov R_inc R_unemp R_lag R_c , nocons;

Iteration 0:   log likelihood = -450872.83  
Iteration 1:   log likelihood = -115395.44  
Iteration 2:   log likelihood = -112234.29  
Iteration 3:   log likelihood = -111873.01  
Iteration 4:   log likelihood = -111872.48  
Iteration 5:   log likelihood = -111872.48  

Logistic regression                             Number of obs     =    650,472
                                                Wald chi2(8)      =  219672.53
Log likelihood = -111872.48                     Prob > chi2       =     0.0000

------------------------------------------------------------------------------
     nc_year |      Coef.   Std. Err.      z    P>|z|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
      Lmajor |   -.807559   .0120894   -66.80   0.000    -.8312537   -.7838642
     R_black |   .0003664   .0010453     0.35   0.726    -.0016823    .0024151
      R_hisp |   .0095177   .0010106     9.42   0.000      .007537    .0114984
       R_pov |  -.0386205   .0042623    -9.06   0.000    -.0469746   -.0302665
       R_inc |  -1.03e-07   1.59e-06    -0.06   0.948    -3.22e-06    3.02e-06
     R_unemp |   .0032477   .0049294     0.66   0.510    -.0064139    .0129092
       R_lag |   8.356694   .0242077   345.21   0.000     8.309248    8.404141
         R_c |  -6.642356   .1133675   -58.59   0.000    -6.864552   -6.420159
------------------------------------------------------------------------------

.         matrix B1 = e(b);

. /*matrices of yr and state fes*/
> matrix yr = (0,0,0,0,0,0,0,0,0,0,0,0);

. matrix st = (0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0);

. #delimit;
delimiter now ;
.         /* Strategic logit with partial observability */
> ml model lf slpo_lf (A_L: nc_year = major , nocons) 
>                                         (A_Rr:  = PercentBlack PercentHisp PercentPoverty MedHsInc PercentUnemploy L.nc_year) 
>                                         (B_Rr:  = public PercentBlack PercentHisp PercentPoverty MedHsInc PercentUnemploy ) , tech(nr) ;

. ml init B1 B2, copy ;

. ml maximize, diff;

initial:       log likelihood = -195689.84
rescale:       log likelihood = -172848.59
rescale eq:    log likelihood =  -122544.1
Iteration 0:   log likelihood =  -122544.1  (not concave)
Iteration 1:   log likelihood = -112376.07  
Iteration 2:   log likelihood = -111019.28  
Iteration 3:   log likelihood = -110733.75  
Iteration 4:   log likelihood = -110713.45  
Iteration 5:   log likelihood = -110710.35  
Iteration 6:   log likelihood = -110710.23  
Iteration 7:   log likelihood = -110710.23  

                                                Number of obs     =    650,472
                                                Wald chi2(1)      =    6058.40
Log likelihood = -110710.23                     Prob > chi2       =     0.0000

---------------------------------------------------------------------------------
        nc_year |      Coef.   Std. Err.      z    P>|z|     [95% Conf. Interval]
----------------+----------------------------------------------------------------
A_L             |
          major |  -1.184275    .015215   -77.84   0.000    -1.214095   -1.154454
----------------+----------------------------------------------------------------
A_Rr            |
   PercentBlack |   .0009892   .0010923     0.91   0.365    -.0011517    .0031301
    PercentHisp |   .0093277   .0010169     9.17   0.000     .0073346    .0113209
 PercentPoverty |    -.03951   .0045082    -8.76   0.000     -.048346    -.030674
       MedHsInc |  -1.46e-06   1.61e-06    -0.91   0.362    -4.61e-06    1.69e-06
PercentUnemploy |  -.0041497   .0051146    -0.81   0.417    -.0141741    .0058748
                |
        nc_year |
            L1. |   8.616734   .2435318    35.38   0.000      8.13942    9.094047
                |
          _cons |  -4.151502   .1163638   -35.68   0.000    -4.379571   -3.923434
----------------+----------------------------------------------------------------
B_Rr            |
         public |    .004223    .034577     0.12   0.903    -.0635465    .0719926
   PercentBlack |  -.0035419   .0008511    -4.16   0.000    -.0052101   -.0018737
    PercentHisp |  -.0067464   .0007839    -8.61   0.000    -.0082828   -.0052101
 PercentPoverty |    -.02165    .003553    -6.09   0.000    -.0286138   -.0146862
       MedHsInc |  -1.10e-06   1.36e-06    -0.81   0.420    -3.76e-06    1.57e-06
PercentUnemploy |   .0401714   .0041652     9.64   0.000     .0320078    .0483351
          _cons |   1.434743   .0995235    14.42   0.000      1.23968    1.629805
---------------------------------------------------------------------------------

.         matrix ll_slpo = e(ll);

. gen sample_m = 1 if e(sample)==1;
(54,206 missing values generated)

.                                                 /* Model fit */
>                         /* i's likelihood */
>                         gen double xb_A_L = ([A_L]_b[major]*major );

.                         gen double xb_A_Rr = ([A_Rr]_b[PercentBlack]*PercentBlack + [A_Rr]_b[PercentHisp]*PercentHisp + [A_Rr]_b[Percent
> Poverty]*PercentPoverty +
>                                                                                 [A_Rr]_b[MedHsInc]*MedHsInc + [A_Rr]_b[PercentUnemploy]*
> PercentUnemploy + [A_Rr]_b[L.nc_year]*L.nc_year + [A_Rr]_b[_cons]);
(54,206 missing values generated)

.                         gen double xb_B_Rr = ([B_Rr]_b[public]*public + [B_Rr]_b[PercentBlack]*PercentBlack + [B_Rr]_b[PercentHisp]*Perc
> entHisp + [B_Rr]_b[PercentPoverty]*PercentPoverty +
>                                                                                 [B_Rr]_b[MedHsInc]*MedHsInc + [B_Rr]_b[PercentUnemploy]*
> PercentUnemploy + [B_Rr]_b[_cons]);
(54,206 missing values generated)

.                         gen double A_L = exp(xb_A_L)/(exp(xb_A_L)+exp(xb_A_Rr));
(54,206 missing values generated)

.                         gen double A_Rr = exp(xb_A_Rr)/(exp(xb_A_L)+exp(xb_A_Rr));
(54,206 missing values generated)

.                         gen double B_Rr = invlogit(xb_B_Rr);
(54,206 missing values generated)

.                                                                         gen double lnfj_slpo = .;
(704,678 missing values generated)

.                         replace lnfj_slpo = ln(A_L + A_Rr*(1-B_Rr)) if nc_year == 0 & sample_m==1;
(591,895 real changes made)

.                         replace lnfj_slpo = ln(A_Rr*B_Rr) if nc_year == 1 & sample_m==1;
(58,577 real changes made)

.                                                 gen pr_san_slpo = A_Rr*B_Rr;
(54,206 missing values generated)

.                         gen pr_det_slpo = A_L;
(54,206 missing values generated)

.                         gen pr_acq_slpo = A_Rr*(1-B_Rr);
(54,206 missing values generated)

.                         gen pred_vio_slpo = 0 if sample_m==1;
(54,206 missing values generated)

.                                 replace pred_vio_slpo =1 if pr_san_slpo>max(pr_det_slpo,pr_acq_slpo) & sample_m==1;
(52,932 real changes made)

.                         gen corr_slpo = 0 if sample_m==1;
(54,206 missing values generated)

.                                 replace corr_slpo = 1 if (pred_vio_slpo==1 & nc_year==1) | (pred_vio_slpo==0 & nc_year==0);
(618,091 real changes made)

.                         gen corr_vio_slpo = .;
(704,678 missing values generated)

.                                 replace corr_vio_slpo = 0 if (pred_vio_slpo==0 & nc_year==1);
(19,013 real changes made)

.                                 replace corr_vio_slpo = 1 if (pred_vio_slpo==1 & nc_year==1);
(39,564 real changes made)

.                         gen incorr_vio_slpo =.;
(704,678 missing values generated)

.                                 replace incorr_vio_slpo = 0 if (pred_vio_slpo==0 & nc_year==0);
(578,527 real changes made)

.                                 replace incorr_vio_slpo = 1 if (pred_vio_slpo==1 & nc_year==0);
(13,368 real changes made)

.                                         sum corr_slpo corr_vio_slpo incorr_vio_slpo;

    Variable |        Obs        Mean    Std. Dev.       Min        Max
-------------+---------------------------------------------------------
   corr_slpo |    650,472    .9502192    .2174919          0          1
corr_vio_s~o |     58,577    .6754187    .4682222          0          1
incorr_vio~o |    591,895    .0225851    .1485767          0          1

.                                                         gen pred_Y3_slpo = A_Rr*(1-B_Rr);
(54,206 missing values generated)

.                         gen Acq_gt_33 = 0 if sample_m==1;
(54,206 missing values generated)

.                                 replace Acq_gt_33 = 1 if pr_acq_slpo>0.30 & sample_m==1;
(1,253 real changes made)

.                         tab Acq_gt_33;

  Acq_gt_33 |      Freq.     Percent        Cum.
------------+-----------------------------------
          0 |    649,219       99.81       99.81
          1 |      1,253        0.19      100.00
------------+-----------------------------------
      Total |    650,472      100.00

.                                                 /* Voung Test*/
>                         svmat ll_l;

.                         svmat ll_slpo;

.                         gen num_V = ll_slpo1 - ll_l1 + ((15/2)*ln(Obs)-(9/2)*ln(Obs));
(704,677 missing values generated)

.                         egen denom1_V_1 = mean((lnfj_slpo/lnfj_l)^2);

.                         egen denom1_V_2 = mean(lnfj_slpo/lnfj_l);

.                         gen denom1_V_3 = denom1_V_2^2;

.                         gen denom2_V = sqrt(denom1_V_1 - denom1_V_3);

.                         gen Voung = num_V/denom2_V;
(704,677 missing values generated)

.                         sum num_V denom2_V Voung;

    Variable |        Obs        Mean    Std. Dev.       Min        Max
-------------+---------------------------------------------------------
       num_V |          1    1189.531           .   1189.531   1189.531
    denom2_V |    704,678    .2774704           0   .2774704   .2774704
       Voung |          1    4287.058           .   4287.058   4287.058

.                                                 /* Clarke test */
>                         signtest lnfj_slpo=lnfj_l;

Sign test

        sign |    observed    expected
-------------+------------------------
    positive |      603207      325236
    negative |       47265      325236
        zero |           0           0
-------------+------------------------
         all |      650472      650472

One-sided tests:
  Ho: median of lnfj_slpo - lnfj_l = 0 vs.
  Ha: median of lnfj_slpo - lnfj_l > 0
      Pr(#positive >= 603207) =
         Binomial(n = 650472, x >= 603207, p = 0.5) =  0.0000

  Ho: median of lnfj_slpo - lnfj_l = 0 vs.
  Ha: median of lnfj_slpo - lnfj_l < 0
      Pr(#negative >= 47265) =
         Binomial(n = 650472, x >= 47265, p = 0.5) =  1.0000

Two-sided test:
  Ho: median of lnfj_slpo - lnfj_l = 0 vs.
  Ha: median of lnfj_slpo - lnfj_l != 0
      Pr(#positive >= 603207 or #negative >= 603207) =
         min(1, 2*Binomial(n = 650472, x >= 603207, p = 0.5)) =  0.0000

.                                                                 /* Summary fit stats */
>                         sum corr_slpo corr_vio_slpo incorr_vio_slpo;

    Variable |        Obs        Mean    Std. Dev.       Min        Max
-------------+---------------------------------------------------------
   corr_slpo |    650,472    .9502192    .2174919          0          1
corr_vio_s~o |     58,577    .6754187    .4682222          0          1
incorr_vio~o |    591,895    .0225851    .1485767          0          1

.                         sum corr_l corr_vio_l incorr_vio_l;

    Variable |        Obs        Mean    Std. Dev.       Min        Max
-------------+---------------------------------------------------------
      corr_l |    650,472    .9502192    .2174919          0          1
  corr_vio_l |     58,577    .6754187    .4682222          0          1
incorr_vio_l |    591,895    .0225851    .1485767          0          1

.                         sum pred_Y3_slpo;

    Variable |        Obs        Mean    Std. Dev.       Min        Max
-------------+---------------------------------------------------------
pred_Y3_slpo |    650,472    .0215689    .0619277   .0006566   .4248438

.                         tab Acq_gt_33;

  Acq_gt_33 |      Freq.     Percent        Cum.
------------+-----------------------------------
          0 |    649,219       99.81       99.81
          1 |      1,253        0.19      100.00
------------+-----------------------------------
      Total |    650,472      100.00

.                 log off;
      name:  <unnamed>
       log:  C:\Users\mdnieman\Desktop\SS Rep\SS_SMPO_epa.log
  log type:  text
 paused on:  29 May 2018, 16:10:44
------------------------------------------------------------------------------------------------------------------------------------------
